Load packages
library("Thermimage")
library("fields")
library("imager")
library("magrittr")
library("spatstat")
Set configurations
## Noise filering
NF = TRUE
## Histgoram Equalization
HE = TRUE
## INPUT_PATH_FOLDER
input = "../Data/"
## OUTPUT_PATH_FOLDER
output = "4_HE_NF/"
## Gassian sigma value for noise filtering
sigma = 0.8
## Situational parameters
OD = 10
RH = 53
AtmosT = 28
Load path images from directory
filenames <- list.files(input, pattern="*.JPG", full.names=FALSE)
filenames
[1] "DJI_0001_R.JPG" "DJI_0003_R.JPG" "DJI_0005_R.JPG" "DJI_0007_R.JPG" "DJI_0009_R.JPG" "DJI_0011_R.JPG"
[7] "DJI_0013_R.JPG" "DJI_0015_R.JPG" "DJI_0017_R.JPG" "DJI_0019_R.JPG" "DJI_0021_R.JPG" "DJI_0023_R.JPG"
[13] "DJI_0025_R.JPG" "DJI_0027_R.JPG" "DJI_0029_R.JPG" "DJI_0031_R.JPG" "DJI_0033_R.JPG" "DJI_0035_R.JPG"
[19] "DJI_0037_R.JPG" "DJI_0039_R.JPG" "DJI_0041_R.JPG" "DJI_0043_R.JPG" "DJI_0045_R.JPG" "DJI_0047_R.JPG"
[25] "DJI_0049_R.JPG" "DJI_0051_R.JPG" "DJI_0053_R.JPG" "DJI_0055_R.JPG" "DJI_0057_R.JPG" "DJI_0059_R.JPG"
[31] "DJI_0061_R.JPG" "DJI_0063_R.JPG" "DJI_0065_R.JPG" "DJI_0067_R.JPG" "DJI_0069_R.JPG" "DJI_0071_R.JPG"
[37] "DJI_0073_R.JPG" "DJI_0075_R.JPG" "DJI_0077_R.JPG" "DJI_0079_R.JPG" "DJI_0081_R.JPG" "DJI_0083_R.JPG"
[43] "DJI_0085_R.JPG" "DJI_0087_R.JPG" "DJI_0089_R.JPG" "DJI_0091_R.JPG" "DJI_0093_R.JPG" "DJI_0095_R.JPG"
[49] "DJI_0097_R.JPG" "DJI_0099_R.JPG" "DJI_0101_R.JPG" "DJI_0103_R.JPG" "DJI_0105_R.JPG" "DJI_0107_R.JPG"
[55] "DJI_0109_R.JPG" "DJI_0111_R.JPG" "DJI_0113_R.JPG" "DJI_0115_R.JPG" "DJI_0117_R.JPG" "DJI_0119_R.JPG"
[61] "DJI_0121_R.JPG" "DJI_0123_R.JPG" "DJI_0125_R.JPG" "DJI_0127_R.JPG" "DJI_0129_R.JPG" "DJI_0131_R.JPG"
[67] "DJI_0133_R.JPG" "DJI_0135_R.JPG" "DJI_0137_R.JPG" "DJI_0139_R.JPG" "DJI_0141_R.JPG" "DJI_0143_R.JPG"
[73] "DJI_0145_R.JPG" "DJI_0147_R.JPG" "DJI_0149_R.JPG" "DJI_0151_R.JPG" "DJI_0153_R.JPG" "DJI_0155_R.JPG"
[79] "DJI_0157_R.JPG" "DJI_0159_R.JPG" "DJI_0161_R.JPG" "DJI_0163_R.JPG" "DJI_0165_R.JPG" "DJI_0167_R.JPG"
[85] "DJI_0169_R.JPG" "DJI_0171_R.JPG" "DJI_0173_R.JPG" "DJI_0175_R.JPG" "DJI_0177_R.JPG" "DJI_0179_R.JPG"
[91] "DJI_0181_R.JPG" "DJI_0183_R.JPG" "DJI_0185_R.JPG" "DJI_0187_R.JPG" "DJI_0189_R.JPG" "DJI_0191_R.JPG"
[97] "DJI_0193_R.JPG" "DJI_0195_R.JPG" "DJI_0197_R.JPG" "DJI_0199_R.JPG" "DJI_0201_R.JPG" "DJI_0203_R.JPG"
[103] "DJI_0205_R.JPG" "DJI_0207_R.JPG" "DJI_0209_R.JPG" "DJI_0211_R.JPG" "DJI_0213_R.JPG" "DJI_0215_R.JPG"
[109] "DJI_0217_R.JPG" "DJI_0219_R.JPG" "DJI_0221_R.JPG" "DJI_0223_R.JPG" "DJI_0225_R.JPG" "DJI_0227_R.JPG"
[115] "DJI_0229_R.JPG" "DJI_0231_R.JPG" "DJI_0233_R.JPG" "DJI_0235_R.JPG" "DJI_0237_R.JPG" "DJI_0239_R.JPG"
[121] "DJI_0241_R.JPG" "DJI_0243_R.JPG" "DJI_0245_R.JPG" "DJI_0247_R.JPG" "DJI_0249_R.JPG" "DJI_0251_R.JPG"
[127] "DJI_0253_R.JPG" "DJI_0255_R.JPG" "DJI_0257_R.JPG" "DJI_0259_R.JPG" "DJI_0261_R.JPG" "DJI_0263_R.JPG"
[133] "DJI_0265_R.JPG" "DJI_0267_R.JPG" "DJI_0269_R.JPG" "DJI_0271_R.JPG" "DJI_0273_R.JPG" "DJI_0275_R.JPG"
[139] "DJI_0277_R.JPG" "DJI_0279_R.JPG" "DJI_0281_R.JPG" "DJI_0283_R.JPG" "DJI_0285_R.JPG" "DJI_0287_R.JPG"
[145] "DJI_0289_R.JPG" "DJI_0291_R.JPG" "DJI_0293_R.JPG" "DJI_0295_R.JPG" "DJI_0297_R.JPG" "DJI_0299_R.JPG"
[151] "DJI_0301_R.JPG" "DJI_0303_R.JPG" "DJI_0305_R.JPG" "DJI_0307_R.JPG" "DJI_0309_R.JPG" "DJI_0311_R.JPG"
[157] "DJI_0313_R.JPG"
Read image function
# Create a function with arguments.
new.read_image <- function(image_file,
input_path="", output_path="Output/",
NF=FALSE,
sigma=0.8,
HE=FALSE,
OD=NULL,
RH=NULL,
AtmosT=NULL) {
## Create image path
image_path = paste(input_path, image_file, sep="")
## Extract meta-tags from thermal image file ##
cams<-flirsettings(image_path, exiftool="installed", camvals="")
## Set variables for calculation of temperature values from raw A/D sensor data
Emissivity<-cams$Info$Emissivity # Image Saved Emissivity - should be ~0.95 or 0.96
ObjectEmissivity<-0.96 # Object Emissivity - should be ~0.95 or 0.96
dateOriginal<-cams$Dates$DateTimeOriginal
dateModif<- cams$Dates$FileModificationDateTime
PlanckR1<- cams$Info$PlanckR1 # Planck R1 constant for camera
PlanckB<- cams$Info$PlanckB # Planck B constant for camera
PlanckF<- cams$Info$PlanckF # Planck F constant for camera
PlanckO<- cams$Info$PlanckO # Planck O constant for camera
PlanckR2<- cams$Info$PlanckR2 # Planck R2 constant for camera
ATA1<- cams$Info$AtmosphericTransAlpha1 # Atmospheric attenuation constant
ATA2<- cams$Info$AtmosphericTransAlpha2 # Atmospheric attenuation constant
ATB1<- cams$Info$AtmosphericTransBeta1 # Atmospheric attenuation constant
ATB2<- cams$Info$AtmosphericTransBeta2 # Atmospheric attenuation constant
ATX<- cams$Info$AtmosphericTransX # Atmospheric attenuation constant
if (is.null(OD)) {
OD<- cams$Info$ObjectDistance # object distance in metres
}
FD<- cams$Info$FocusDistance # focus distance in metres
ReflT<- cams$Info$ReflectedApparentTemperature # Reflected apparent temperature
if (is.null(AtmosT)) {
AtmosT<- cams$Info$AtmosphericTemperature # Atmospheric temperature
}
IRWinT<- cams$Info$IRWindowTemperature # IR Window Temperature
IRWinTran<- cams$Info$IRWindowTransmission # IR Window transparency
if (is.null(RH)) {
RH<- cams$Info$RelativeHumidity # Relative Humidity
}
h<- cams$Info$RawThermalImageHeight # sensor height (i.e. image height)
w<- cams$Info$RawThermalImageWidth # sensor width (i.e. image width)
## Rotate image
temp_image = rotate270.matrix(readflirJPG(image_path))
## Convert raw to temp data
temp_image = raw2temp(temp_image, ObjectEmissivity, OD, ReflT, AtmosT, IRWinT, IRWinTran, RH,
PlanckR1, PlanckB, PlanckF, PlanckO, PlanckR2)
## Noise filtering with gauss filter
if (NF == TRUE) {
temp_image = as.im(temp_image)
temp_image = blur(temp_image, sigma, bleed=FALSE)
temp_image = as.matrix(temp_image)
}
# hist(temp_image, main=paste("Histogram ", image_file))
## Save temperature image
saveRDS(temp_image, file = paste(output_path, tools::file_path_sans_ext(image_file), ".RData", sep=""))
## Return either min or max or every temperature value
if (HE == TRUE) {
temp_image
} else {
min <- min(temp_image)
max <- max(temp_image)
values <- c(min, max)
}
}
Load images
glob_temp_range <- unlist(lapply(filenames, new.read_image,
input_path=input,
output_path=output,
NF=NF,
sigma=sigma,
HE=HE,
OD=OD,
RH=RH,
AtmosT=AtmosT))
Calculate global min and max values
min_temp <- min(glob_temp_range)
max_temp <- max(glob_temp_range)
Calculate ECDF
if (HE == TRUE) {
f <- ecdf(glob_temp_range)
hist(glob_temp_range, main=paste("Global Histogram"))
hist(f(glob_temp_range), main=paste("Global Equalized Histogram"))
} else {
f <- NULL
}
Normalize images
new.normalize <- function(image_file, output_path="4_HE_NF/",
min, max,
HE=FALSE, f_ecdf=NULL) {
## Create image path
image_path = paste(output_path, tools::file_path_sans_ext(image_file), ".RData", sep="")
temp_image <- readRDS(image_path)
# unlink(image_path)
## Normalize data
temp_image <- mirror.matrix(temp_image)
if (HE == TRUE) {
equalized <- f_ecdf(temp_image)
equalized <- as.cimg(equalized, dim=dim(temp_image))
plot(equalized, rescale=FALSE, main=paste("Histogram equalized", image_file))
save.image(equalized, file=paste(output_path, image_file, sep=""), quality = 1)
} else {
temp_image <- (temp_image - min) / (max - min)
# hist(x=temp_image, main=paste("Histogram of normalized values", image_path), xlab='normalized value')
## Plot and save temperature image
temp_image = as.cimg(temp_image, dim=dim(temp_image))
save.image(temp_image, file=paste(output_path, image_file, sep=""), quality = 1)
plot(temp_image, rescale=FALSE)
}
}
Plot normalized images
temp = lapply(filenames, new.normalize,
output_path=output,
min=min_temp, max=max_temp,
HE=HE,
f_ecdf=f)